#include <stdlib.h>
#include <stdio.h>
#include <gmp.h>
#include <math.h>
#include <pthread.h>

#define N_ITERATES 24
#define PRECISION 35000000

mpf_t a, a_, b, b_,t, t_, p, p_, aux, pi;

pthread_t calc_a, calc_b, calc_t, calc_p, calc_pi;

void *proximo_a(void *v){
	mpf_add(a, a_, b_);	// a(n+1)=an+bn 
	mpf_div_ui(a, a, 2); // a(n+1)=a(n+1)/2 
}

void *proximo_b(void *v){
	mpf_mul(b, a_, b_);	// b(n+1)=an*bn 
	mpf_sqrt(b, b);	// b(n+1)=raiz(b(n+1)) 
}

void *proximo_t(void *v){
	mpf_sub(t, a_, a);	// t(n+1)=an-a(n+1) 
	mpf_pow_ui(t, t, 2);	// t(n+1)=t(n+1)²
	mpf_mul(t, p_, t);	// t(n+1)=pn*t(n+1)
	mpf_sub(t, t_, t);	// t(n+1)=tn-t(n+1) 
}

void *proximo_p(void *v){
	mpf_mul_ui(p, p, 2); // p(n+1)=2*pn 
}

void *proximo_pi(void *v){
	mpf_add(pi, a, b);	// pi=a(n+1)+b(n+1) 
	mpf_pow_ui(pi, pi, 2);	// pi=pi² 
	mpf_mul_ui(aux, t, 4);	// aux=4*t(n+1) 
	mpf_div(pi, pi, aux);	// pi=pi/aux */
}

int main(int argc, const char *argv[]) {

	int i;
	time_t start, end;

	mpf_set_default_prec(PRECISION);
	gmp_printf("precisao em bits: %d\n", (int)mpf_get_default_prec());
	mpf_inits(a_, t_, b_, p_, aux, pi, NULL);
	mpf_init_set_d(a, 1.0);
	mpf_init_set_d(t, 0.25);
	mpf_init_set_d(b, (1/sqrt(2)));
	mpf_init_set_ui(p, 1);

	start = time(0);

	for (i = 0; i < N_ITERATES; i++) {
		mpf_set(a_, a);   // a_=an 
		mpf_set(t_, t);	// t_=tn 
		mpf_set(b_, b);	// b_=bn 
		mpf_set(p_, p);	// p_=pn 
		pthread_create(&calc_a, NULL, proximo_a, NULL);	//thread para calcular o proximo valor de a
		pthread_create(&calc_b, NULL, proximo_b, NULL);	//thread para calcular o proximo valor de b
		pthread_join(calc_a, NULL); //para calcular t, o proximo valor de a precisa ser calculado primeiro
		pthread_create(&calc_t, NULL, proximo_t, NULL); //thread para calcular o proximo valor de t
		pthread_create(&calc_p, NULL, proximo_p, NULL); //thread para calcular o valor de p
		pthread_join(calc_t, NULL); 
		pthread_join(calc_b, NULL); //para calcular pi, o proximo valor de b precisa ser calculado antes
		pthread_create(&calc_pi, NULL, proximo_pi, NULL); //thread para calcular o valor de pi
		pthread_join(calc_p, NULL);	
		pthread_join(calc_pi, NULL);
		//nesse ponto todas as threads terminaram de executar
		gmp_printf("%da interacao: %.10000000Ff\n\n", i+1, pi);
	}

	end = time(0);

	printf("Tempo de calculo: %d segundos\n", end-start);

	mpf_clears(a, a_, b, b_, t, t_, p, p_, aux, pi, NULL);

	return 0;
}
